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Abstract 



We consider the liquid-gas phase boundary in a binary fluid mixture near 
its critical end point. Using general scaling arguments we show that the 
diameter of the liquid-gas coexistence curve exhibits singular behaviour as 
the critical end point is approached. This prediction is tested by means of 
extensive Monte-Carlo simulations of a symmetrical Lennard-Jones binary 
mixture within the grand canonical ensemble. The simulation results show 
clear evidence for the proposed singularity, as well as confirming a previously 
predicted singularity in the coexistence chemical potential [Fisher and Upton, 
Phys. Rev. Lett. 65, 2402 (1990)]. The results suggest that the observed 
singularities, particularly that in the coexistence diameter, should also be 
detectable experimentally. 
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I. INTRODUCTION 



A critical end point occurs when a line of second order phase transitions intersects and is 
truncated by a first order phase boundary, beyond which a new non- critical phase is formed. 
Critical end points are common features in the phase diagrams of a many multicomponent 
systems such as binary fluid mixtures, binary alloys and liquid crystals. They also occur 
under certain circumstances in pure systems having additional internal degrees of freedom 
such as superfiuids, as well as certain ferromagnets and ferroelectrics. 

Perhaps the simplest system to exhibit a critical end point is a binary fluid mixture. 
Here the phase diagram is spanned by three thermodynamic fields: the temperature T, a 
chemical potential /i, and an ordering field h coupling to the concentration difference of the 
two particle species. A number of different phase diagram topologies are possible for binary 
fluids depending on the microscopic interaction parameters, as systematically classified by 
Konynenburg and Scott The particular case on which we shall focus in this work is 
of type II in their classification scheme, and is depicted schematically in figure |I| for the 
subspace h = 0. Within this subspace, two fluid phases (3 and 7, each rich in one of the two 
particle species, coexist with one another. By tuning T and /i, however, one finds a critical 
'A' line, Tc{fi), where both phases merge into a single disordered jSj phase. The point at 
which the A-line Tc(/i) intersects the first order line of liquid-gas transitions fia{T) marks 
the critical end point (Te,/Xe)- For T < Te, the phase boundary fia(T) constitutes a triple 
line along which the fluid phases (3 and 7 coexist with the gas phase a, while for T > Te, 
fia{T) defines the region where the /37 and a phases coexist. Precisely at the critical end 
point the critical mixture of jS and 7 phases coexists with the gas phase. Since the gas phase 
does not participate in the criticality, it is commonly referred to as a "spectator" phase. 

Interest in critical end points has recently been aroused following theoretical work by 
Fisher and coworkers [0. On the basis of phenomenological scaling and thermodynamic 
arguments, these authors argued that the nonanalytic behaviour at the critical end point 
engenders a free energy like singularity in the first order (spectator) phase boundary /Xo-(T). 
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The nature of this singularity, specifically its amplitude ratio, was related to the universal 
features of the critical line. Additionally, new universal amplitude ratios were proposed 
for the noncritical surface tensions near the critical end point 0. These predictions were 
subsequently corroborated by analytical calculations on extended spherical models as 
well as a Landau theory study of a model ferroelectric 0. To date, however, experimental 
results pertaining to critical end points are scarce and those that do exist have focused on the 
interfacial properties and the surface amplitude ratios To the best of our knowledge 

there have been no reported attempts to study the bulk coexistence properties and the 
matter of the predicted singularity in the spectator phase boundary. There is a similar 
dearth of simulation work on the subject, and we know of no detailed numerical studies of 
critical end point behaviour, either in lattice or continuum models. 

In this paper we address the issue of the nature of the spectator phase boundary near a 
critical end point by means of computer simulation of a continuum binary fluid model 0. 
Our paper is organised as follows. In section |I| we briefly review the principal features of 
the scaling arguments of Fisher and coworkers, and show that in addition to the previously 
predicted singularity in fi^{T), they also imply singular behaviour in the diameter of the 
liquid-gas coexistence curve at the critical end point. In section |ITI| we detail extensive 
Monte-Carlo simulations of a symmetrical (Lennard- Jones) binary fluid mixture. We map 
the liquid-vapour coexistence curve and the A-line of the model by applying finite-size scaling 
analyses to the probability distribution functions of appropriate observables. The results 
provide remarkably clear signatures of divergences in the appropriate temperature derivatives 
of the coexistence diameter and the phase boundary chemical potential, thus corroborating 
the theoretical predictions. Finally in section |^ we detail our conclusions. 



II. THE LIQUID-GAS COEXISTENCE CURVE 
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A. The coexistence chemical potential 

In this subsection, we briefly review the scahng arguments of Fisher and Barbosa p| 
concerning the singularity in the spectator phase boundary at the critical end point. In so 
doing, we shall continue to employ the language of the binary fluid mixture, although the 
arguments themselves are not restricted to this case. 

Within the grand canonical framework, liquid-gas coexistence is prescribed by the equal- 
ity of the Gibbs free energy G = — ksT In Z in the respective phases i.e. 

Gg{fi^{T),T,h) =Gi{fi^{T),T,h) (2.1) 

Since the gas spectator phase is necessarily noncritical, its free energy is analytic at the end 
point and can thus be expanded as 

Ggifi, T, h)=Ge + GfAfi + Gft + G^h + GfAfi^ + ■■■ (2.2) 

where t = T — Tg, and Afi = ^ — /^e, with A/i < 0. 

The liquid phase on the other hand is critical and therefore contains both an analytic 
(background) and a singular contribution to the free energy 

G'K/i,T,/i) = Go(/i,T,/i) - |rp-"e;±(/i|r|-^) (2.3) 

where Go is the analytic part, while G^{y) is a universal scaling function which is a function of 
the relevant scaling fields r(T, fi, h) = T — Tc{fi) and h{T, fi,h) ^ h that measure deviations 
from the A-line. Q^{y) must satisfy matching conditions as y ^ ±oo, and the quantities a 
and A are respectively the specific heat and gap exponents associated with the A-line. 
To linear order the scaling fields may be expanded as 

t(T, fi, h) = t + aih + a2Afi (2.4a) 
h{T,ix,h) = h + bit + b2Afx (2.4b) 



where the and 6j are non universal 'field mixing' parameters [pi?| , P | some of which have 
geometrical significance within the phase diagram. In particular, since r = specifies Tc{fi), 
one has at the critical end point 



representing the gradient of the A-hne at the critical end point. 



The critical free energy, equation |2.3|, can also be expanded in A/i, t and h. Recalling 



equation |2.2|, invoking equation and solving for yu(T) then yields 



/i^(T) - iio{T) ^ -X±|t|2-" - Y±\tf\h\ - ^Z±\t\-^h\ (2.6) 

valid as T — >^ Te±, — > 0. Here X±, Y±, Z± are critical amplitudes, and a, (3, 7 are the usual 
critical exponents characterising the A-line. Corrections to scaling have been neglected and 
the background term has the expansion 

/io(T) = fJ^e + git + g2h ^ (2.7) 

If we restrict out attention to the coexistence surface h = and assume a > 0, equa- 
tion implies a specific heat-like divergence in the curvature of the spectator phase bound- 
ary 

^^-^±1^'^ (2.8) 
where the amplitude ratio X^/X^ is expected to be universal 0. 

B. The coexistence diameter 

Let us now examine the behaviour of the coexistence particle number density in the 
vicinity of the critical end point. Specifically, we shall focus on the coexistence diameter in 
the region h = 0, defined by 

P,iT) ^ ^[pgifx^iT)) + piif^^iT))]. (2.9) 



The particle density is obtainable from the Gibbs potential as 

1 /dG 
'V \dfl 



1 dG\ 

P = -T7{ — ] (2-10) 

T,h 



so that 



1 (dGMT)^ 

^^^^ = -V [ d-, , 



(2.11) 



with 



a(/i^(T),T) = [Gg{fi,{T),T) + Gi{fi,{T),T)]/2. 



(2.12) 



Appeahng to eqs. |2.2| and p.3| , and noting that /5 = 2 — a — A, then yields for the singular 
behaviour 



Pd{T) = U±\T{T,fi^{T))f + V±\T{T,fi^{T))\'-" 
+ terms analytic at Tg 

as r — ^ 0. Here the non-universal critical amplitudes take the form 

14 = a2(2-a)6;±(0) 

with Q±{z) = dQ±/dz. 



Now, along the liquid-gas coexistence curve, one has from equation |2.4a 



|r(T,/x^(T))| = \t + a2{i^a{T) -/ie)|. 



Recalhng equations |2]^ and (and setting /i = 0), one then finds 



|r(r,/x^(r))| = |t|[l + a2^i + 0(|t|i-")]. 



Thus to leading order one can write 



Pd{T) = lJ±\tf + VM-^ 

-|- terms analytic at Tg, 



(2.13) 



(2.14a) 
(2.14b) 



(2.15) 



(2.16) 



(2.17) 



where f/± = (1 + a2giYU± and V± = (1 + a2giY °'V±. We note that this expression is of the 
same form as the singularity in the overall density on the A— line [1T0 |. 
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A special case of equation |2.17|, relevant to the present work, is that for a symmetrical 



fluid having energetic invariance under h ^ —h. In this case one finds on symmetry grounds 



that the field mixing parameters 61 = 62 = and hence from equation p. 14a 



p^(T) = V±\t\'"" + terms analytic at T^, (2.18) 
which implies a divergent diameter derivative 

^^V^\tr (2.19) 

where Kt = (1 — a)V±. Since this divergence occurs in the first derivative of the observable 
Pd{T), it is in principle more readily visible than that in the second derivative of ficr{T), cf. 
eq. p.8| . As we shall now show, however, clear signatures of both divergences are readily 
demonstrable by Monte-Carlo simulation. 



III. MONTE CARLO SIMULATIONS 

A. Model and algorithmic considerations 

The simulations described here were performed for a symmetrical binary fluid model 
using a Metropolis algorithm within the grand canonical ensemble (GCE) ||r^. The fluid 
is assumed to be contained in volume V = with periodic boundary conditions. The 
grand-canonical partition function takes the form 

00 00 N /■ „ N 

= E E n I / d^i] e['^^-*«^^»+'^(^-^^)l (3.1) 

Ni=0N2=0i=l ^ ^ 

where $ = J2i<j (Pi^ij) is the total configurational energy, /i is the chemical potential, and h 
is the ordering field (all in units of ksT). N = Ni + N2 is the total number of particles of 
types 1 and 2. 

The interaction potential between particles i and j was assigned the familiar Lennard- 
Jones (LJ) form 

0(r.,) = 4e„„[(a/r,,)^2 _ ^cr/r,,)% (3.2) 
7 



where a is a parameter which serves to set the interaction range, while emn measures the 
well-depth for interactions between particles of types m and n. In common with most 
other simulations of Lennard- Jones systems, the potential was truncated at a cutoff radius 
Tc = 2.5a (irrespective of the species of the interacting particles) and left unshifted. To 
simplify the task of locating interacting particles, a cell decomposition scheme was employed 
in which the total system volume was partitioned into cubic cells of side Tc- Interaction 
emanating from a specific particle then only extend as far as the 26 neighbouring cells. 

An Ising like symmetry was imposed on the model by choosing en = €22 = e > 0. This 
choice endows the system with energetic invariance under h —h, thereby ensuring that 
the critical end point lies in the surface h = 0. We shall accordingly restrict our attention 
henceforth to this regime. A further parameter ei2 = 6 was used to control interactions 
between unlike particles. The phase diagram of the model in the surface h = is then 
uniquely parameterised by the ratio 6/e. Physically, the role of this ratio in determining 
the form of the phase diagram can be understood as follows. For 6/e ^ 1, the energy 
penalty associated with contacts between unlike particles is small and hence there is little 
incentive for phase separation unless the temperature is very low and the density very high. 
One therefore expects a phase diagram having a critical end point temperature <^ Tj:^ 
and density pe ^ Pc^, where Tj:^ and are the liquid-gas critical temperature and density 
respectively. Choosing a smaller value of 6/e, however, moves the end point towards the 
liquid-gas critical point, into which it merges for sufficiently small 6/e, forming a tricritical 



point [|T3|. Empirically we find that for 6/e = 0.6, there is a tricritical point, while for 
6/e = 0.75 there is a critical end point having pe ~ 2.3p^^. The location of the critical end 
point for larger values of 6/e could not be reliably determined since its density is above 
that at which the grand canonical particle insertion algorithm (see below) is operable. In 
this work, all simulations were performed with 6/e = 0.7, which yields critical end point 
parameters Tg ^ 0.93T^'^, pe ~ 1.75p^^. This temperature is sufficiently small compared to 
Tj.^ that critical density fluctuations don't obscure the end point behaviour, while at the 
same time pe is not so large as to unduely hinder particle insertions. 
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In order to sample ergodically the phase space of the model, two sorts of Monte-Carlo 
action are necessary. The first is a particle transfer step in which one attempts either to 
insert a particle at a randomly chosen position, or alternatively, to delete a randomly chosen 
existing particle. Candidate particles for insertion are assigned a species type (1 or 2) with 
equal probability. The second sort of action is an identity swop, in which one chooses an 
existing particle at random and attempts to change its identity (1 ^ 2 or 2 — 1). Combined 
use of these operations samples the grand canonical ensemble in which the particle density 
p = N/V, energy density u = ^/V and number difference density m = {Ni — N2)/V all 
fluctuate. 

In accordance with convention, we shall employ dimensionless units to express our data: 

p = pa^, fh = ma'^, u = ua^ (3.3) 

f = ksT/e (3.4) 

We also note for future reference that our algorithm utilises not the true chemical potential 
/i, but an effective chemical potential p to which the true chemical potential is related by 

p = p + Po-\n [n/L^) (3.5) 

where po is the chemical potential in the non- interacting (ideal gas) limit. It is this effective 
value that features in the results that follows. 



B. Method and Results 

In the course of the simulations, three systems sizes of volume V = (T.Scr)^, V = (lOu)^ 
and V = (12. 5o")'^ were studied, corresponding to average particle numbers ~ 250, ~ 
600 and ^ 1200 respectively at the critical end point (whose location we discuss below). 
Following equilibration, runs comprising up to 6 x 10^ MCS |]T^ were performed and the 
joint distribution pi(p, m, -u) was sampled approximately every 10"^ MCS and accumulated 
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in the form of a histogram. For each L, simulations were carried out at several (typically 5) 
temperatures along the liquid-gas coexistence curve. 

In order to locate the liquid-gas coexistence curve, the finite-size form of the density 
distribution Pl{p) was studied. Precisely at coexistence, Pl{p) is (to within corrections ex- 
ponentially small in L) double peaked with equal weight in both peaks For a given 



simulation temperature, this 'equal weight' criterion can be used to determine the coexis- 
tence chemical potential to high accuracy. Close to the liquid-gas critical point i.e. when 
T ^ Tj:^, this task is straightforward; one simple tunes the chemical potential at constant 
T until the measured form of Pl{p) is double peaked. By contrast, however, in the strongly 
first order regime (T <C Tj:^) , this approach cannot be used to obtain the coexistence form of 
Pl{p) due to the large free energy barrier separating the coexisting phases. This barrier leads 
to metastability effects and prohibitively long correlation times. To circumvent this diffi- 
culty one must employ advanced sampling schemes. One such scheme is the multicanonical 
preweighting method ||T6[ which encourages the simulation to sample the interfacial config- 
urations of intrinsically low probability. This is achieved by incorporating a suitably chosen 
weight function in the Monte-Carlo update probabilities. The weights are subsequently 
'folded out' from the sampled distribution to yield the correct Boltzmann distributed quan- 
tities. Use of this method permits the direct measurement of the order parameter probability 
distribution at first order phase transitions, even when the distribution spans many decades 
of probability. 

To fully capitalise on the data gathered in the simulations, the histogram reweighting 



technique [|T§] was employed. This technique rests on the observation that histograms accu- 
mulated at one set of model parameters can be reweighted to yield estimates of histograms 
appropriate to another set of not-too-distant parameters. Results from individual simulation 
runs at different parameters can also be combined in a systematic fashion to provide infor- 
mation over larger regions of the phase diagram. When used in tandem with multicanonical 
preweighting, the histogram reweighting technique constitutes a powerful tool for mapping 
the coexistence curve properties of continuum fluid models fl^ . 
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Using the multicanonical preweighting technique, simulations of the hquid-gas coexis- 
tence curve of the binary fluid were performed. To begin with, the region near the liquid-gas 
critical point was studied since there the free energy barrier between the coexisting phases 
is small and no preweighting function is need. Thereafter, the temperature was reduced in 
a stepwise fashion and the coexistence form of Pl{p) was collected. For each successively 
lower temperature studied, histogram reweighting was used to obtain a suitable preweight- 
ing function by extrapolating from the coexistence histograms previously obtained at higher 



temperatures. Further details concerning this strategy have been given elsewhere Fig- 
ure IKa) shows the resulting coexistence forms of Pl{,p) for the system of size L = lOa, at 
a selection of temperatures in the range 0.9T^^ <T < T^a. Also shown in figure |(b) is the 
same data expressed on a log scale. The corresponding estimates for the coexistence chem- 
ical potential and the coexistence diameter are plotted in figures ^ and ^ respectively. The 
gas and liquid densities required for the diameter calculation were obtained as the average 
densities of the respective peaks of Pl{p)- It is interesting to note from figure H(b), that 
for the lowest temperature studied (T ^ 0.9Tc), the ratio of minimum to maximum in the 
distribution is approximately 10^^. A free energy barrier of this magnitude would, of course, 
represent an insurmountable obstacle to phase space evolution, were it not for the use of 
multicanonical preweighting. 

It is also instructive to examine the coexistence behaviour of the number difference 
order parameter distribution ^^(m), for temperatures above and below Tg, as shown in 
figure 1^. Well below Tg, this distribution is three-peaked, with one narrow peak centered on 
m = and two broader peaks centered on positive and negative values of m. The peak at 
m = corresponds to the disordered gas phase, while the degenerate peaks at positive and 
negative values of m represent the ordered A-rich and i?-rich liquid phases respectively. As 
one approaches the critical end point, however, the liquid peaks become broader and overlap 
with the gas phase peak, ultimately leaving just one single peak at m = for T > Tg. 

To determine the locus of the A-line and the critical end point in which it formally 
terminates, finite-size scaling (FSS) methods were employed. On the A-line, criticality is 
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signalled by the asymptotic scale invariance of pL{rri). A useful dimensionless measure of 
the form of this distribution is its fourth order cumulant ratio Ul = 1 — 3(m'^) / (m^)^ ||T^. If 
one plots Ul{T) for a given constant /2, one expects the curves for each L to intersect at the 
critical de-mixing temperature. For binary fluids with short ranged interactions, the critical 
behaviour on the A-line and the critical end point |]2D| is expected to be Ising like. One thus 



expects that asymptotically (for sufficiently large L) there will be a unique crossing point 
at a universal value U^. Ising model simulations |2l|] give the estimate = 0.470(3). 



We have measured Ul{T) for several value of the chemical potential fi, thus enabling an 
estimate of the locus of the A-line. For each chosen value of fl, an initial rough estimate of the 
critical temperature was obtained from a number of short simulations, in which the form of 
Pii'fn) was observed visually. Longer runs at this temperature were then carried out for each 
system size in order to facilitate a more precise determination of Tc. A representative plot of 
Ul as a function of T and L is shown in figure |^ for the case fl = —2.95. Each of the curves 
shown represent the histogram extrapolation of data obtained from a single simulation at 
the temperature T = 1.005. Evidently there is no unique cumulant intersection point for 
the three system sizes, rather the crossing point for the L = 7.5a and L = lOcr system sizes 
occurs at larger T and smaller Ul than for the L = 10a and L = 12.5cr system sizes. It is 
reasonable to assume, however, that for sufficiently large L, the cumulant intersection point 
would converge on the universal value, = 0.470(3), but that for the system sizes studied 
here this limit is not attained due to corrections to finite-size scaling and (to a lesser extent) 
statistical error. In view of this, we adopt as our estimate of T'c(/i), that temperature at 
which the cumulant intersection occurs for the two largest system sizes. 

The critical end point was estimated in a similar manner, by studying the cumulant 
ratio Ul for Pl{^) as a function of T and L along the liquid branch of the coexistence curve. 
In principle, the form of Pl(^^) corresponding to the liquid phase can be extracted from 
the joint coexistence distribution ^^(p, m,£t), since this embodies contributions from both 
phases. Assuming there is a clear distinction between the gas and liquid phases (ie. when 
T ^ Tj:^), one can isolate the liquid branch component simply by disregarding all those 
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histogram entries having p < pd- However, it turns out to be difficult to obtain sufficient 
statistics for ^^(m) in this way. This is because the preweighted simulation uniformly 
samples all densities between the gas and liquid phases. When the gas and liquid densities 
are well separated, the simulation spends only a relatively small fraction of the time sampling 
liquid-like configurations and accordingly the statistical quality of the histogram for ^^(m) 
is low. To circumvent this problem, it is expedient to perform a simulation that samples only 
the liquid phase, while neverthless remaining at coexistence. This is achievable by setting 
T and jl to their coexistence values and performing a simulation without preweighting the 
density distribution. If a starting configuration having a liquid-like density is chosen, the 
simulation will then remain in the liquid phase by virtue of the large free energy barrier 
separating it from the gas phase. The results of implementing this procedure are shown 
in figure |^, in which f/^ is plotted as a function of temperature along the liquid branch of 
the coexistence curve in the neighbourhood of the critical end point. The curves exhibit an 
intersection at = 0.958(2), with a cumulant value Ul = 0.47 for the two largest system 
sizes and Ul = 0.455 for the smallest two system sizes. The associated estimate of the 
critical end point chemical potential and density are pe = —3. 017(3), pe = 0.587(5). 

In the interests of completeness we have also obtained the location of the liquid-gas 
critical point. This was achieved via a FSS study of the density-like ordering operator. 
Details of the approach used have been extensively described elsewhere [^ , p^ and so will 
not be repeated here. We merely quote the results: f^.^ = 1.024(2), p[f = 0.327(2), = 
-2.767(4). 

In figure ^ we present the estimated phase diagram of our system. Plotted are the 
coexistence densities for the liquid and gas branches, the position of the liquid-gas critical 
point, the locus of the A-line and the position of the critical end point. One observes that 
close to the liquid-gas critical point, strong finite-size effects occur in the peak densities of 
Pl{p)- These have also been observed and discussed in the context of a previous study of the 
pure LJ fluid |T^. More striking, however, is the appearance of a pronounced 'kink' in the 
liquid branch density close to the critical end point. This is a consequence of the singular 

13 



behaviour on the critical hne. As one would expect, the gas branch displays no such kink 
due to the analyticity of Gg{fi,T) at Tg. 

To probe more closely the behaviour of the coexistence density, we plot in figure ^(a) the 
diameter derivative —dpd{T)/dT, for the three system sizes studied. The data exhibit a clear 
peak close to Tg, that narrows and grows with increasing system size. Very similar behaviour 
is also observed in the curvature of the spectator phase boundary —d?jla/dT'^, see figure 
11(b). These peaks constitute, we believe, the finite-size-rounded forms of the divergences 
eqs. p.l9| and On the basis of finite-size scaling theory [|l9l, the peaks are expected to 
grow in height like L"/'^, with u the correlation length exponent. Unfortunately, it is not 
generally feasible to extract estimates of a/v in this way (even for simulations of lattice 
Ising models), because to do so requires the ability to measure the analytic background, for 
which the present system sizes are much too small. Nevertheless, the correspondence of the 
peak position with the independently estimated value of Te, as well as the narrowing and 
growth of the peak with increasing L constitutes strong evidence supporting the existence 
of the predicted singularities. 



IV. DISCUSSION 

In summary, we have employed grand canonical Monte- Carlo simulations in conjunction 
with multicanonical preweighting and histogram reweighting to study the first order phase 
boundary near the critical end point of a continuum binary fluid model. The results provide 
clear evidence for the existence of singularities in the phase boundary chemical potential and 
the coexistence curve diameter. They thus constitute strong corroboration of the scaling 
arguments on whose basis the singularities are predicted. 

Although no experimental observations of singularities at critical end points have yet 
been reported, we believe that given an appropriately chosen system their presence should 
be relatively easy to detect. In this respect, a binary fluid model might be a good candidate 
system for study, since generally speaking rather precise liquid-gas coexistence curve data 
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can be obtained. In this case the role of the chemical potential in equation p.6| is replaced 
by the pressure, but otherwise the form of the singularity remains unaltered. Of course in a 
real system, the lack of finite-size rounding should render the singularities more conspicuous 
than in a simulation. Moreover, since real binary fluids do not generally possess a special 
symmetry between the two fluid components, the chemical potential and temperature feature 
in the scaling field /i, changing the diameter singularity into a much stronger l^l'' 

singularity (cf. section In future work we also intend to study the consequences of this 
symmetry breaking, via simulations of an asymmetrical binary fluid model. 
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FIG. 1. Schematic phase diagram of a binary fluid in the coexistence surface /i = 0. The broken 
line /Xo-(T) is the first order hquid-gas phase boundary between the fluid and gas phase a. The full 
line is the critical line of second order transitions T\{[i) separating the demixed phases /?+7, from 
the mixed phase /37. The two lines intersect at the critical end point. 
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FIG. 2. (a). Estimates of the coexistence density distributions for the L = 10a systems size, 
for a range of subcritical temperatures, obtained as described in the text. The hnes are merely 
guides to the eye. Statistical errors do not exceed the symbol sizes, (b). The same data expressed 
on a log scale. 
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FIG. 3. Estimates of the the number difference order parameter distribution pL{rh) for the 
L = lOcr system size at three temperatures spanning Tg. 
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FIG. 4. The temperature dependence of the coexistence chemical potential jl obtained from 
the histogram reweighting of the distributions shown in figure |2[ Statistical errors do not exceed 
the symbol sizes. 
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FIG. 5. The temperature dependence of the coexistence curve diameter pd = {pg + pi)/2^ 
obtained from histogram reweighting of the distributions shown in figure |2[ Statistical errors 
not exceed the symbol sizes. 
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FIG. 6. The measured cumulant ratio Ul{T) for the three system sizes studied. The Unes are 
the results of histogram reweighting of single simulations performed at T = 1.005, p, = —2.95 
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FIG. 7. The cumulant ratio Ul for the distribution pi,{rh), obtained according to the procedure 
described in the text. The hues are the results of histogram reweighting of hquid-phase simulations 
performed at liquid-gas coexistence. 
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FIG. 8. The average peak densities corresponding to the coexistence form of pl{p) for the three 
systems sizes studied, plotted as a function of temperature. Also shown is the estimated locus of 
the A-line (filled circles) and the liquid-gas critical point (cross). Statistical errors do not exceed 
the symbol sizes. 
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FIG. 9. (a) The numerical temperature derivative of the measured coexistence diameter 
—dpd/dT in the vicinity of the critical end point temperature, (b) The measured curvature of 
the phase boundary — defter /dT'^ in the vicinity of the critical end point temperature. In both cases 
data are shown for the three system sizes studied. 
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